Quantum dynamics of repulsively bound atom pairs in the Bose-Hubbard model 



00 

o 
o 



in 

(N 



5h 

43 



o3 



C 

O 

(N 
> 

O 
in 

o 

o 

-i— > 



T3 

c 

o 
o 



X 



Li Wang, Yajiang Hao, and Shu CherQ 
Beijing National Laboratory for Condensed Matter Physics, 
Institute of Physics, Chinese Academy of Sciences, Beijing 100080, China 
(Dated: March 25, 2008) 

We investigate the quantum dynamics of repulsively bound atom pairs in an optical lattice de- 
scribed by the periodic Bose-Hubbard model both analytically and numerically. In the strongly 
repulsive limit, we analytically study the dynamical problem by the perturbation method with the 
hopping terms treated as a perturbation. For a finite-size system, we numerically solve the dynamic 
problem in the whole regime of interaction by the exact diagonalization method. Our results show 
that the initially prepared atom pairs are dynamically stable and the dissociation of atom pairs is 
greatly suppressed when the strength of the on-site interaction is much greater than the tunneling 
amplitude, i.e., the strongly repulsive interaction induces a self-localization phenomenon of the atom 
pairs. 

PACS numbers: 03.75.Kk, 03.75.Lm, 32.80.Pj 



I. INTRODUCTION 

The ultracold bosonic atoms in optical lattices have 
opened a new window to study many-particle quantum 
physics in a uniquely controlled manner [3, S El • Vari- 
ous schemes have been proposed to realize a wide range 
of many-body Hamiltonians by loading ultracold atoms 
into a certain optical lattice @, H, [H • The advances in 
the manipulation of atom-atom interactions by Feshbach 
resonances 7] allowed experimental study of many-body 
systems accessible to the full regime of interactions,even 
to the very strongly interacting Tonks-gas limit [H, 
Recently a lot of exciting experiments in optical lattices 
have been implemented, including the superfluid-Mott- 
insulator transition 0, |a|, n on-linear self-trapping in a 
periodical optical lattice [lOfland repulsively bound atom 
pairs in an optical lattice The basic physics of the 
ultracold atomic systems in optical lattice is captured by 
the Bose-Hubbard model (BHM) 0, i], which incorpo- 
rates the competition between the interaction and hop- 
ping energy and has been successfully applied to interpret 
the superfluid-Mott-insulator transition. As a fundamen- 
tal model, the BHM has been widely applied to study the 
quantum phase transitions and dynamic problems in the 
optical lattices. 

Very recently, Winkler et al. [ll[ have studied the dy- 
namical evolution of the initially prepared bosonic atom 
pairs in an optical lattice. Their experimental results 
indicate that the atom pairs with strongly repulsive on- 
site interaction U exhibits a remarkably longer lifetime 
than the system with weakly repulsive interaction [IT1 |. 
At first glance, this result is counter-intuition because 
one may expect the repulsive interaction to separate the 
particles, instead to bind them together. The experimen- 
tal result has stimulated theoretical investigations on the 
dynamics of repulsively bound pairs [TH, [l]| . In Ref . [HI , 



the theoretical understanding of the stable pair relies on 
the analytical solution of a two-particle problem by solv- 
ing two particle Lippmann-Schwinger scattering equa- 
tion on the lattice corresponding to the Bose-Hubbard 
Hamiltonianpll [l2| . Obviously, this method is only lim- 
ited to a single repulsively bound pair and is not capable 
to extend to deal with many-particle dynamic problem. 

Motivated by the experimental progress in this 
paper we study the quantum dynamics of the repulsively 
pair states in the BHM both analytically and numerically. 
In the strongly repulsive limit, we develop an analytical 
method to deal with the dynamical problem based on the 
perturbation expansion of the hopping terms, whereas 
we can numerically solve the dynamic problem in the 
whole regime of interaction for a finite-size system which 
could be diagonalized exactly by the exact diagonaliza- 
tion method. The Bose-Hubbard Hamiltonian (BHH) 
reads HQ 
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where bl(bj) is the creation (annihilation) operator of 
bosons on site i (j), hi — S|S, counts the number of 
bosons on the ith site, and denotes summation over 
nearest neighbors. The parameter J denotes the hop- 
ping matrix element to neighboring sites, and U cx a s 
represents the on-site interaction due to s-wave scatter- 
ing. For an actual optical lattice, J and U are related to 
the depth of the optical lattice Vo which is determined by 
the intensity of the laser. The lattice constant a is half 
of the wave length of the laser A/2 Q. In this article, 
we focus on the dynamical evolution of the repulsively 
bound atom pairs in the periodic Bose-Hubbard model 
with U > 0, i.e. repulsive on-site interaction. In the fol- 
lowing calculation, we will ignore all possible dissipations 
in the system, such as the loss of atoms by three-body 
collision. 

The paper is organized as follows. In Sec. II, we first 
review a general scheme to deal with dynamical evolution 
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and present the exact result for the two-site problem. In 
Sec. Ill, we develop a perturbation method to study the 
dynamic evolution of the initially prepared state of atom 
pairs which works in the large U limit. In Sec. IV, we 
study the dynamical problem for a finite-size system by 
using the exact diagonalization method and compare the 
analytical results with the exact numerical results. 



II. GENERAL SCHEME 

The Bose-Hubbard model has been investigated by a 
variety of theoretical and numerical methods under dif- 
ferent cases [H, EE US E3- Most of the theoretical in- 
vestigations concern the ground state properties, whereas 
the quantum dynamic problems are hardly dealt with and 
most of works are limited to the double-well (two-site) 
problems [H,E1- 

Given an initial state \ipo) , the evolution state at time 
t can be formally represented as 

m))=e- iSt \ip ). (2) 

If we know all the eigenstates of the BHH which are the 
solutions of schrodinger equation H \ 4> n ) = E n \<j> n ) , we 
can get 

\m)=J2 e ' iBnt ^^Ho). (3) 

n 

Now it is straightforward to obtain the probability of 
finding the given initial state at time t 

p (t) = \(m\m)\ 2 

= EE|a„| 2 |a m | 2 e^^"^, (4) 

n— 1 m— 1 

where a n = {<j) n |^o) and D = {N + L - l)\/[N\ (L - 1)!] 
is the basis dimension of the TV-particle system in a lat- 
tice with size L. If the energy spectrum and its corre- 
sponding eigenstates <j) n are known, then the dynamical 
problem is exactly solved in principle. It is obvious that 
the dimension D increases very quickly as the particle 
number N and lattice site L increase. Therefore it is not 
practical to investigate the dynamical problem of a large 
system in this way, though one can solve it numerically 
for a finite-size system by exact diagonalization method. 

It is instructive to study the two-site system before 
presenting the analytical results based on perturbation 
theory and the numerical results for the systems with 
larger size. The dynamical evolution of an atom pair 
in the toy two-site model can be completely solved an- 
alytically. The results of this toy model will shed light 
on studies of larger systems. According to the scheme 
discussed above, it is quite easy to get the analytical ex- 
pression of Po(t) for a pair of atoms initially in a two-site 




FIG. 1: (Color online) Time evolution of Po(t) versus t in 
units of J for the two-site-two-boson case. A comparison be- 
tween the exact and the approximate solution is made. The 
solid lines denote the exact solutions given in Sec. II. The 
dotted lines denote the approximate solutions by perturbation 
theory. From top to bottom the lines correspond to U / J = 60, 
30 and 10, respectively. 

lattice: 

P (t) = | ai | 4 + |a 2 | 4 + |a 3 | 4 

+2 \ ai \ 2 \a 2 \ 2 cos (Ei - E 2 )t 
+2 |a 2 | 2 |a 3 | 2 cos (E 2 - E 3 )t 
+2\a 3 \ 2 M 2 cos (E 3 - E x )t, (5) 

where a x = 1/V2, a 2 = l/V2 + a 2 and a 3 = 1/V 2 + P 2 
with the parameters a = —4y/2J/ (U - VU 2 + 16 J 2 ) 
and (3 = -4V2J/ (U + VU 2 + 16 J 2 ). The eigenener- 
gies are given by E x = U, E 2 = (U - VU 2 + 16 J 2 ) /2, 

and E 3 = (U + VU 2 + 16J 2 ) /2. It is quite easy to show 
the time evolution of the initial atom pair by the formula 
([5]). In Fig. 1, we plot the time evolution of Po(t) for 
U/J = 60, 30 and 10. The results show that the larger 
the repulsive U is, the probability of finding the initial 
atom pair is larger. When U/J—+ oo, we observe that 
a 2 -> 0, a 3 -> 1/V2, E\ — E 2 w E 3 — E 2 as U and 
E 3 — E\ w 4J 2 /U. In this limit, the terms including of 
a 2 in Eq. |(SJ), which is proportional to y cos(Z7f), merely 
give a very tiny contribution to Pq (t) accompanying with 
a very high-frequency oscillation, whereas the dominant 
part is given by 

Po(t)^l + \co S (4p). (6) 

It is obvious that P (t) » 1 if Jt < U/J. Therefore at 
the large U limit, the atom pair is dynamically stable 
within a very large time scale. 
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III. DYNAMICS IN THE LARGE U LIMIT 



In the strongly repulsive interaction limit of U / J 1, 
it is convenient to represent the BHH as 

H = H[j + Hj, 

where Hjj and Hj are the on-site interaction part and 
tunnelling part respectively. We may solve the following 
dynamical evolution equation 

\m) = e-^-^ l^o) (7) 

approximately by treating Hj as a perturbation. In order 
to do that, we formally represent \ip(t)) as 



-iHut 



U(t) IVo) 



(8) 



with U(t) defined as U(t) = e iHvt e -Wjt-iH„t _ It is eagy 
to check that U(t) fulfills the following motion of equation 



d_ 

df 



U{t) = -iHj(t)U(t) 



with 



Hj(t) = e lHut Hje- lHut . 
Therefore the formal solution of U (t) is given by 



U(t) 



Texp 



o 

t f t 

dti ■ ■ ■ dt n T 
i) Jo 



Hj{tx)---Hj{t r 



where T is the time-ordering operator. By using the 
above formula, now we can calculate (ip(0) \ip{t)) per- 
turbatively. Up to the second order, we get 

wo) \m) 



Ho) - 

(H e- lRut f dt, I' dt 2 Hj{ti)Hj(t 2 ) |Vo) • (9) 
Jo Jo 

From the definition of \ip(t)) , we see that it should be 
normalized. However, in the scheme of perturbation cal- 
culation, the normalized condition is not automatically 
fulfilled. In order to overcome this problem and get a 

well-defined Po(t), we need use a normalized ipit)^ which 

is defined as 



W)) 



vW) \m) 



(10) 



Therefore in the scheme of perturbation expansion, the 
probability Po(t) is defined by 



p (t) = (v(o) m 



(ii) 
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FIG. 2: (Color online) Time evolution of Po(t) for the six-site- 
six-boson case. The comparison between the solution given 
by perturbation theory and the exact solution from numerical 
computation is shown. The analytical results are denoted by 
the dotted lines and the numerical result by solid lines. The 
top two lines correspond to U / J — 60 and the other two lines 
correspond to U / J — 30. The time is in units of J. 



For a given initial state |?/>o) represented in a basis of 
Fock states 



|V>o) = \ni,ri2, ■ 



,n L ) 



(12) 



N U |Vo) 



with rii — or 2, it is easy to get Hjj \tpo) - „ 
For convenience, we first calculate the case with |-i/>o 
|2, 0,2,0..., 2,0). By using the formula ©, we get 



wo) m)) 



-in 



Therefore we have 



2NJH 2NJ 2 



jut 



U 



Po{t) 



A 



1 



U 2 



'AN J 2 87V 2 J 4 



(13) 



V u 2 

AN 2 JH 2 8N 2 J 4 t 



U 2 



IP 



IP 

sin Ut 



(1-cos Ut) + 
(14) 



where A — (ip(t) IV'(i)) is the normalized constant. We 
note that, in the scheme of the second-order perturba- 
tion theory, formula (fT3|) has the same form for different 
initial configurations such as \ipo) = |2, 0, 2, 0, 0, 2, 0, 0, 0) 
and \i/j ) = |2,0,2,0,2,0,0,0,0). Actually, the formula 
(fT5)l is valid for all the initial configurations where the 
atoms pairs are not neighboring to each others. For the 
other kind of initial configurations, it is straightforward 
to calculate -PoW by formula (fTTj) and 

In order to check the validity of the perturbation the- 
ory, we plot the Po(t) obtained by the perturbation the- 
ory for the two-site problem in Fig. [T]and make a rough 
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FIG. 3: Time evolution of Po(t) for the 50-particle system 
with initial state |2, 0, 2, 0, 2, 0, 2, 0) versus t in units of J 
with the corresponding on-site interaction U/J = 200. 



comparison between the perturbation result and the ex- 
act result given in Eq. (|5|). Specially, for the simplest 
two-site case with \ijj } = |2, 0) or |0, 2), we have 



Po(t) 



A- 1 
8JH 



4J 2 
U 2 



8J 4 



(1 - cos Ut) + 



AJH 



4+2 



IP 



sin Ut 



IV. NUMERICAL RESULT 

The exact diagonalization (ED) method is very pow- 
erful from the viewpoint that it can give out all informa- 
tion of the system exactly. To solve the matrix eigenvalue 
problem of the Hamiltonian, it is convenient to work in 
a basis of Fock states 



| a) = |{ni,n 2 , -,n L } a ) 



(17) 



with rii = 0, 1,2, ••• and the index a — 1, ...,£) labels 
the different compositions of occupation numbers for a 



fixed total particle number J2i 



N. Through the 



exact diagonalization algorithm we compute exactly all 
the eigenvalues E n and the corresponding eigenstates <p n 
of the Hamiltonian (Q]). Then we can calculate the dy- 
namic evolution of Po(t) by using Eq. To give a 
concrete example and compare with the result obtained 
by the perturbation theory in Sec. Ill, we first consider 
an initial state \ipo) = |2, 0, 2, 0, 2, 0). Our result is shown 
in Fig. O We can see that at large U limit, the pertur- 
bation theory works quite well. 

We note that what was measured in the experiment 
fill ] is the number of remaining pairs after a variable hold 
time. To understand the experimental result qualita- 
tively, we calculate the time-dependent normalized num- 
ber of the atom pairs which can be represented as 



N p (t) = f- pia)Kit)l2 



M 



(18) 



(15) where c a (t) is given by 



where the normalization constant A is given by 



A = l + 8 J 4 



2 2t t 2 
— r (1 - cos Ut) sin Ut H 

c/ 4 v ' u 3 u 2 



(16) 



As is shown in Fig. [TJ the larger the U, the better the 
exact and the approximate solution fit. 

It is natural to use the formula (fTTj) to deal with the 
many-particle system. From (|13p , we see that the pertur- 
bation theory should work well as long as NJ 2 t/U <C 1. 
In order to compare with the exact result by numerical 
exact diagonalization method in next section, we apply 
(fT4")l to calculate the time evolution of Po (t) for a finite- 
size system with an initial state l^o) = |2, 0, 2, 0, 2, 0) . 
Our result is shown in Fig. [2] We can see that at large 
U limit, the perturbation theory works quite well. Our 
results show that ip(t) decays slower for a larger U. For a 
larger system with N = 50 and the initial state given by 
|2, 0,2,0, ...,2,0,2,0), the time evolution of P (t) shows 
similar behaviors as that of the 2- and 6-particle systems. 
As shown in Fig. O Po(t) decays slowly accompanying 
with a high-frequency but very narrow oscillation. These 
results show that, when U is large enough, tp(t) mainly 
stays in the initial state within a large time scale. The 
larger the U is, the more obvious the case. 



«(*) = £;• 



— iE n t 



(19) 



N p (a) denotes the number of atom pairs in \a) and M = 
N/2 is the total number of initial atom pairs. Using the 
formula (|18[) . we can calculate the normalized atom pair 
number for different on-site interaction U and different 
initial state \ipo). 

For simplicity, let us consider a one-dimensional (ID) 
optical lattice with L = 9 and M — 3, i.e. three boson 
atom pairs loaded into an optical lattice of nine sites. 
The time dependent normalized pair numbers for vari- 
ous interaction strengths are plotted in Fig. Q] with the 
initial state given by |2, 2, 2, 0, 0, 0, 0, 0, 0). Evidently, the 
atom pairs in the case of U = 100.0 stay longer than that 
in the case of U = 1.0. It is clear that the large on-site 
interaction U effectively suppresses the dynamic dissocia- 
tion of the initial pairs, which is in qualitative agreement 
with the experiment reported in Ref. [TTI ] . 

Next we study how the dynamics of the repulsively 
bound pairs is affected by the choice of the initial state. 
For the case with three bound pairs loaded into an op- 
tical lattice with nine sites, there are C| = 84 different 
configurations. Taken the translational symmetry of the 
system into account, there are in fact only seven different 
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FIG. 4: (Color online) Time evolution of normalized atom 
pair number in a nine-site optical lattice for U/J = 1.0, 10.0 
and 100.0. The time is in units of J and the initial state is 
chosen as [2,2,2,0,0,0,0,0,0). 



configurations. For different initial state \ipo)i we calcu- 
late the normalized atom pair number and plot them in 
Fig. [5] According to Fig. [5j although for different ini- 
tial state |^o) the corresponding dynamics of the system 
is similar to each other, there are still some minor dif- 
ferences depending on whether the initial pairs are close 
together or not. For initial states |2, 0, 2, 0, 2, 0, 0, 0, 0), 
|2,0,2,0,0,2,0,0,0), and (2,0,0,2,0,0,2,0,0) in which 
each pair is apart from the others, the time evolutions 
are almost the same. This is consistent with our result 
by perturbation theory. According to the second-order 
perturbation theory, the time evolution of Po(t) for these 
initial states are same. Similarly, the time evolution for 
initial states (2,2,0,2,0,0,0,0,0), (2,2,0,0,2,0,0,0,0), 
and |2, 2, 0, 0, 0, 2, 0, 0, 0) are almost the same for differ- 
ent on-site interaction U as is shown in Fig. [5l We 
notice that the initial states with bound pairs repelling 
each other are dynamically stabler than the states with 
all atom pairs initially occupying neighboring sites. We 
also investigate the case in which there are three boson 
atom pairs loaded into a two-dimensional 3x3 optical 
lattice. The results are similar to the one-dimensional 
case as shown in Fig. 4 and Fig. 5. 

From above discussions, we could understand the sta- 
ble atom-pair state due to the strongly repulsive interac- 
tion as a dynamically stable state. The physical origin 
is essentially that the energy conservation does not allow 
the particle to tunnel into the neighboring sites to save 
the repulsive interaction unless the obtained tunneling 
energy is comparable to the loss of repulsive energy [2(| • 
From theoretical point of view, a similar dynamically sta- 
ble state is expected to appear for the case in which the 
initial state composed of three or more atoms. In order 
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FIG. 5: (Color online) Time evolution of the normalized atom 
pair number for different on-site interaction U and different 
initial states \ipo). From top to bottom the on-site interaction 
U/J= 1.0, 10.0, and 100.0, respectively. In each graph, there 
are seven lines corresponding to seven different initial states. 
The initial states are denoted by different lines shown in the 
table at the bottom. The time is in units of J. 



to clarify this point, we investigate an example case in 
which six atoms were initially loaded into only one site 
of a ID optical lattice with nine sites. The initial state 
is chosen as |0, 0, 0, 0, 6, 0, 0, 0, 0). Given the initial state, 
we then determine the dynamic evolution of the initial 
state on the whole lattice, i.e., the spatial distribution 
of the normalized atom number on each site at different 
time point for different on-site interaction. The time- 
dependent normalized atom number on the j-th site can 
be calculated by 



nj(t) 



D 

E 

a=l 



N 



(20) 



where n" = 



| a). The spatial distributions provide 



us useful information for the expansion dynamics of a 
condensate loaded in the optical lattice. At an incremen- 
tal sequence of time points, we calculate the normalized 
atom number at each site. The results are shown in Fig. 
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FIG. 6: Time evolution of the normalized atom number on 
each site of the optical lattice. We choose the initial state \ipo) 
as |0, 0, 0, 0, 6, 0, 0, 0, 0). The graph on the top corresponds to 
U/J=1.0 and the graph on the bottom to U/J=W.O. The 
time is in units of J. 
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V. SUMMARY 

In this paper, we theoretically investigate the dynam- 
ical stability of atom pairs induced by the strongly re- 
pulsive interaction in an optical lattice both analytically 
and numerically. Firstly, in the large U limit, we cal- 
culate the time evolution of an initial system composed 
of atom pairs in a scheme of the perturbation theory. 
The analytical results show that the initial state of atom 
pairs is dynamically stabler for larger repulsive interac- 
tion. Then by the exact diagonaliztion method, we nu- 
merically study the stability of atom pairs in a finite-size 
system with few particles by calculating the time evolu- 
tion of the normalized atom pair number. Our results 
show that the initial state of atom pairs are dynamically 
stable and the dissociation of atom pairs is greatly sup- 
pressed when the strength of the on-site interaction is 
large. We also compare our numerical result with the re- 
sult based on the perturbation theory and find out that 
they fit quite well in the large U limit. Our results imply 
that the experimentally observed repulsively bound pairs 
can be understood as a dynamically stable state. 
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